!
! Copyright (C) 2000-2013 A. Marini and the YAMBO team
!              https://code.google.com/p/rocinante.org
!
! This file is distributed under the terms of the GNU
! General Public License. You can redistribute it and/or
! modify it under the terms of the GNU General Public
! License as published by the Free Software Foundation;
! either version 2, or (at your option) any later version.
!
! This program is distributed in the hope that it will
! be useful, but WITHOUT ANY WARRANTY; without even the
! implied warranty of MERCHANTABILITY or FITNESS FOR A
! PARTICULAR PURPOSE.  See the GNU General Public License
! for more details.
!
! You should have received a copy of the GNU General Public
! License along with this program; if not, write to the Free
! Software Foundation, Inc., 59 Temple Place - Suite 330,Boston,
! MA 02111-1307, USA or visit http://www.gnu.org/copyleft/gpl.txt.
!
subroutine io_elemental(ID,VAR,VAR_SZ,CHECK,WARN,OP,I0,DB_I0,I1,DB_I1,&
&                       DB_R0,R0,R1,DB_R1,CH0,DB_CH0,CH1,L0,DB_L0,UNIT,MENU,&
&                       DESCRIPTION,DO_NOT_DUMP)
 !
 use pars,    ONLY:SP,schlen,lchlen
#if defined _NETCDF_IO
 use pars,    ONLY:nf90_SP
 use netcdf
#endif
 use drivers, ONLY:list_dbs
 use com,     ONLY:msg
 use stderr,  ONLY:string_pack,intc,real2ch
 use IO_m,    ONLY:write_is_on,io_unit,netcdf_call,io_netcdf_var,ver_is_gt_or_eq,&
&                  netcdf_dim,io_rec_pos,read_is_on,io_mode,VERIFY,&
&                  io_status,DUMP,REP,LOG,io_com,io_file,io_netcdf_support,&
&                  io_code_version,io_serial_number,io_code_revision
 implicit none
 integer               :: ID
 integer,parameter     :: MAX_DB_ELEMENTS=100
 character(*),optional :: VAR,OP(:)
 integer,     optional :: VAR_SZ,MENU
                               ! MENU = 1 (on) 0 (off)
 integer,     optional :: I0,DB_I0,I1(:),DB_I1(:)
 real(SP),    optional :: R0,DB_R0,R1(:),DB_R1(:),UNIT
 logical,     optional :: L0,DB_L0,CHECK,WARN,DO_NOT_DUMP
 character(*),optional :: CH0,DB_CH0,CH1(:)
 character(lchlen),optional :: DESCRIPTION
 !
 ! Work Space
 !
 character(lchlen) :: ch
 character(1)      :: msg_where
 integer :: i_db(MAX_DB_ELEMENTS),i_err
 real(SP):: r_db(MAX_DB_ELEMENTS),local_unit
 character(schlen) :: db_ch(MAX_DB_ELEMENTS)
 logical :: CHECK_,DUMP_
 !
 DUMP_=io_mode(ID)==DUMP
 if (present(DO_NOT_DUMP)) then
   DUMP_=.not.DO_NOT_DUMP
 endif
 !
 if  (present(CHECK)) then
   CHECK_=CHECK
 else
   CHECK_=.FALSE.
 endif
 !
 local_unit=1._SP
 if (present(UNIT)) local_unit=UNIT
 msg_where=""
 if (io_com(ID)==REP) msg_where="r"
 if (io_com(ID)==LOG.or.list_dbs) msg_where="s"
 !
 if (present(VAR).and.present(VAR_SZ)) then
   io_rec_pos(ID)=1
   if (VAR_SZ>0) then
     if (write_is_on(ID)) then
       if (io_netcdf_support(ID)) then
#if defined _NETCDF_IO
         call netcdf_call(nf90_redef(io_unit(ID)))
         if (present(CH0).or.present(CH1)) then
           call netcdf_call(nf90_def_var(io_unit(ID),VAR,&
&                           nf90_char,(/netcdf_dim(ID,schlen),&
&                           netcdf_dim(ID,VAR_SZ)/),io_netcdf_var(ID)))
         else
           if (ver_is_gt_or_eq(ID,revision=723)) then 
             call netcdf_call(nf90_def_var(io_unit(ID),VAR,&
&                           nf90_SP,netcdf_dim(ID,VAR_SZ),io_netcdf_var(ID)))
           else
             call netcdf_call(nf90_def_var(io_unit(ID),VAR,&
&                           nf90_real,netcdf_dim(ID,VAR_SZ),io_netcdf_var(ID)))
           endif
         endif
         call netcdf_call(nf90_enddef(io_unit(ID))) 
#endif
       endif
       ch=string_pack("[WR",trim(io_file(ID)),"]",&
&                     repeat('-',max( 60-len_trim(io_file(ID))-4,1) ))
     else if (read_is_on(ID)) then
       ch=string_pack("[RD",trim(io_file(ID)),"]",&
&                     repeat('-',max( 60-len_trim(io_file(ID))-4,1) ))
       if (io_netcdf_support(ID)) then
#if defined _NETCDF_IO
         call netcdf_call(nf90_inq_varid(io_unit(ID),VAR,io_netcdf_var(ID)))
#endif
       endif
     endif
     if (.not.present(MENU)) call msg(msg_where,trim(ch))
     if (present(MENU)) then
       if (MENU==1) call msg(msg_where,trim(ch))
     endif
   else
     write (ch,'(a,i6.6,2a,3(a,i2.2),a,i4.4,a)')  '- S/N ',io_serial_number(ID),' ',&
&          repeat('-',26),'- v.',io_code_version(ID,1),'.',io_code_version(ID,2),&
&          '.',io_code_version(ID,3),' r.',io_code_revision(ID),' -'
     if (.not.present(MENU)) call msg(msg_where,trim(ch))
     if (present(MENU)) then
       if (MENU==1) call msg(msg_where,trim(ch))
     endif
   endif
   return
 endif
 !
 if (present(I0)) then
   if (write_is_on(ID)) then
     if (io_netcdf_support(ID)) then
#if defined _NETCDF_IO
       call netcdf_call(nf90_put_var(io_unit(ID),io_netcdf_var(ID),&
&                       I0,(/io_rec_pos(ID)/)))
#endif
     else
       write(io_unit(ID)) I0
     endif
     if (present(DB_I0)) DB_I0=I0
     if (present(VAR)) call msg(msg_where,VAR,I0)
   else if (read_is_on(ID)) then
     if (io_netcdf_support(ID)) then
#if defined _NETCDF_IO
       call netcdf_call(nf90_get_var(io_unit(ID),io_netcdf_var(ID),&
&                       i_db(1),(/io_rec_pos(ID)/)))
#endif
     else
       read(io_unit(ID)) i_db(1)
     endif
     if (DUMP_)          I0=i_db(1)
     if (present(DB_I0)) DB_I0=i_db(1)
     ch=""
     if (io_mode(ID)==VERIFY.and.present(OP)) then
       i_err=i_verify((/I0/),(/i_db(1)/))
       if (CHECK_) then
         if (io_status(ID)==0) io_status(ID)=i_err
         if (i_err/=0.and.present(VAR)) ch=string_pack("*ERR* ",VAR)
       else if (present(WARN)) then
         if (WARN.and.i_err/=0.and.present(VAR)) ch=string_pack("*WRN* ",VAR)
       endif
     endif
     if (len_trim(ch)>0) then
       call msg(msg_where,trim(ch),i_db(1))
     else if (present(VAR)) then
       call msg(msg_where,VAR,i_db(1))
     endif
   endif
   !
   ! DESCRIPTION
   !
   if (present(DESCRIPTION)) then
     DESCRIPTION=VAR//intc(I0)
   endif
   !
   io_rec_pos(ID)=io_rec_pos(ID)+1
   !
 endif
 !
 if (present(I1)) then
   if (write_is_on(ID)) then
     if (io_netcdf_support(ID)) then
#if defined _NETCDF_IO
       call netcdf_call(nf90_put_var(io_unit(ID),io_netcdf_var(ID),&
&                       I1,(/io_rec_pos(ID)/)))
#endif
     else
       write(io_unit(ID)) I1
     endif
     if (present(DB_I1)) DB_I1=I1
     if (present(VAR)) call msg(msg_where,VAR,I1)
   else if (read_is_on(ID)) then
     if (io_netcdf_support(ID)) then
#if defined _NETCDF_IO
       call netcdf_call(nf90_get_var(io_unit(ID),io_netcdf_var(ID),&
&                       i_db(:size(I1)),(/io_rec_pos(ID)/)))
#endif
     else
       read(io_unit(ID)) i_db(:size(I1))
     endif
     if (DUMP_)          I1=i_db(:size(I1))
     if (present(DB_I1)) DB_I1=i_db(:size(I1))
     ch=""
     if (io_mode(ID)==VERIFY.and.present(OP)) then
       i_err=i_verify(I1,i_db(:size(I1)))
       if (CHECK_) then
         if (io_status(ID)==0) io_status(ID)=i_err
         if (i_err/=0.and.present(VAR)) ch=string_pack("*ERR* ",VAR)
       else if (present(WARN)) then
         if (WARN.and.i_err/=0.and.present(VAR)) ch=string_pack("*WRN* ",VAR)
       endif
     endif
     if (len_trim(ch)>0) then
       call msg(msg_where,trim(ch),i_db(:size(I1)))
     else if (present(VAR)) then
       call msg(msg_where,VAR,i_db(:size(I1)))
     endif
   endif
   !
   ! DESCRIPTION
   !
   if (present(DESCRIPTION)) then
     DESCRIPTION=VAR//trim(intc(I1(1)))
     if (size(I1)==2) DESCRIPTION=VAR//trim(intc(I1(1)))//' - '//trim(intc(I1(2)))
   endif
   !
   io_rec_pos(ID)=io_rec_pos(ID)+size(I1)
 endif
 !
 if (present(R0)) then
   if (write_is_on(ID)) then
     if (io_netcdf_support(ID)) then
#if defined _NETCDF_IO
       call netcdf_call(nf90_put_var(io_unit(ID),io_netcdf_var(ID),&
&                       R0,(/io_rec_pos(ID)/)))
#endif
     else
       write(io_unit(ID)) R0
     endif
     if (present(DB_R0)) DB_R0=R0
     if (present(VAR)) call msg(msg_where,VAR,R0*local_unit)
   else if (read_is_on(ID)) then
     if (io_netcdf_support(ID)) then
#if defined _NETCDF_IO
       call netcdf_call(nf90_get_var(io_unit(ID),io_netcdf_var(ID),&
&                       r_db(1),(/io_rec_pos(ID)/)))
#endif
     else
       read(io_unit(ID)) r_db(1)
     endif
     if (DUMP_)          R0=r_db(1)
     if (present(DB_R0)) DB_R0=r_db(1)
     ch=""
     if (io_mode(ID)==VERIFY.and.present(OP)) then
       i_err=r_verify((/R0/),(/r_db(1)/))
       if (CHECK_) then
         if (io_status(ID)==0) io_status(ID)=i_err
         if (i_err/=0.and.present(VAR)) ch=string_pack("*ERR* ",VAR)
       else if (present(WARN)) then
         if (WARN.and.i_err/=0.and.present(VAR)) ch=string_pack("*WRN* ",VAR)
       endif
     endif
     if (len_trim(ch)>0) then
       call msg(msg_where,trim(ch),r_db(1)*local_unit)
     else if (present(VAR)) then
       call msg(msg_where,VAR,r_db(1)*local_unit)
     endif
   endif
   !
   ! DESCRIPTION
   !
   if (present(DESCRIPTION)) then
     DESCRIPTION=VAR//trim(real2ch(R0*local_unit))
   endif
   !
   io_rec_pos(ID)=io_rec_pos(ID)+1
   !
 endif
 !
 if (present(R1)) then
   if (write_is_on(ID)) then
     if (io_netcdf_support(ID)) then
#if defined _NETCDF_IO
       call netcdf_call(nf90_put_var(io_unit(ID),io_netcdf_var(ID),&
&                       R1,(/io_rec_pos(ID)/)))
#endif
     else
       write(io_unit(ID)) R1
     endif
     if (present(DB_R1)) DB_R1=R1
     if (present(VAR)) call msg(msg_where,VAR,R1*local_unit)
   else if (read_is_on(ID)) then
     if (io_netcdf_support(ID)) then
#if defined _NETCDF_IO
       call netcdf_call(nf90_get_var(io_unit(ID),io_netcdf_var(ID),&
&                       r_db(:size(R1)),(/io_rec_pos(ID)/)))
#endif
     else
       read(io_unit(ID)) r_db(:size(R1))
     endif
     if (DUMP_)          R1=r_db(:size(R1))
     if (present(DB_R1)) DB_R1=r_db(:size(R1))
     ch=""
     if (io_mode(ID)==VERIFY.and.present(OP)) then
       i_err=r_verify(R1,r_db(:size(R1)))
       if (CHECK_) then
         if (io_status(ID)==0) io_status(ID)=i_err
         if (i_err/=0.and.present(VAR)) ch=string_pack("*ERR* ",VAR)
       else if (present(WARN)) then
         if (WARN.and.i_err/=0.and.present(VAR)) ch=string_pack("*WRN* ",VAR)
       endif
     endif
     if (len_trim(ch)>0) then
       call msg(msg_where,trim(ch),r_db(:size(R1))*local_unit)
     else if (present(VAR)) then
       call msg(msg_where,VAR,r_db(:size(R1))*local_unit)
     endif
   endif
   !
   ! DESCRIPTION
   !
   if (present(DESCRIPTION)) then
     DESCRIPTION=VAR//trim(real2ch(R1(1)*local_unit))
     if (size(R1)==2) DESCRIPTION=VAR//trim(real2ch(R1(1)*local_unit))//&
&                                 ' - '//trim(real2ch(R1(2)*local_unit))
     if (size(R1)==3) DESCRIPTION=VAR//trim(real2ch(R1(1)*local_unit))//&
&                                 ' '//trim(real2ch(R1(2)*local_unit))//&
&                                 ' '//trim(real2ch(R1(3)*local_unit))
   endif
   !
   io_rec_pos(ID)=io_rec_pos(ID)+size(R1)
 endif
 !
 if (present(L0)) then
   i_db(1)=0
   if (L0) i_db(1)=1
   if (write_is_on(ID)) then
     if (io_netcdf_support(ID)) then
#if defined _NETCDF_IO
       call netcdf_call(nf90_put_var(io_unit(ID),io_netcdf_var(ID),&
&                       i_db(1),(/io_rec_pos(ID)/)))
#endif
     else
       write(io_unit(ID)) i_db(1)
     endif
     if (present(VAR)) call msg(msg_where,VAR,i_db(1)==1)
     if (present(DB_L0)) DB_L0=i_db(1)==1
   else if (read_is_on(ID)) then
     if (io_netcdf_support(ID)) then
#if defined _NETCDF_IO
       call netcdf_call(nf90_get_var(io_unit(ID),io_netcdf_var(ID),&
&                       i_db(2),(/io_rec_pos(ID)/)))
#endif
     else
       read(io_unit(ID)) i_db(2)
     endif
     ch=""
     if (present(DB_L0)) DB_L0=i_db(2)==1
     if (DUMP_)          L0=i_db(2)==1
     if (io_mode(ID)==VERIFY.and.present(OP)) then
       i_err=i_verify((/i_db(1)/),(/i_db(2)/))
       if (CHECK_) then
         if (io_status(ID)==0) io_status(ID)=i_err
         if (i_err/=0.and.present(VAR)) ch=string_pack("*ERR* ",VAR)
       else if (present(WARN)) then
         if (WARN.and.i_err/=0.and.present(VAR)) ch=string_pack("*WRN* ",VAR)
       endif
     endif
     if (len_trim(ch)>0) then
       call msg(msg_where,trim(ch),i_db(2)==1)
     else if (present(VAR)) then
       call msg(msg_where,VAR,i_db(2)==1)
     endif
   endif
   !
   ! DESCRIPTION
   !
   if (present(DESCRIPTION)) then
     DESCRIPTION=VAR//' no'
     if (L0) DESCRIPTION=VAR//' yes'
   endif
   !
   io_rec_pos(ID)=io_rec_pos(ID)+1
 endif
 !
 if (present(CH0)) then
   if (write_is_on(ID)) then
     db_ch(1)=CH0
     if (io_netcdf_support(ID)) then
#if defined _NETCDF_IO
       call netcdf_call(nf90_put_var(io_unit(ID),io_netcdf_var(ID),&
&                       db_ch(:1),(/io_rec_pos(ID)/)))
#endif
     else
       write(io_unit(ID)) db_ch(1)
     endif
     if (present(DB_CH0))   DB_CH0=CH0
     if (present(VAR)) call msg(msg_where,VAR,trim(db_ch(1)))
   else if (read_is_on(ID)) then
     if (io_netcdf_support(ID)) then
#if defined _NETCDF_IO
       call netcdf_call(nf90_get_var(io_unit(ID),io_netcdf_var(ID),&
&                       db_ch(:1),(/io_rec_pos(ID)/)))
#endif
     else
       read(io_unit(ID)) db_ch(1)
     endif
     ch=""
     if (DUMP_)             CH0=db_ch(1)
     if (present(DB_CH0))   DB_CH0=db_ch(1)
     if (io_mode(ID)==VERIFY) then
       i_err=0
       if (trim(db_ch(1))/=CH0) i_err=-1
       if (CHECK_) then
         if (io_status(ID)==0) io_status(ID)=i_err
         if (i_err/=0.and.present(VAR)) ch=string_pack("*ERR* ",VAR)
       else if (present(WARN)) then
         if (WARN.and.i_err/=0.and.present(VAR)) ch=string_pack("*WRN* ",VAR)
       endif
     endif
     if (len_trim(ch)>0) then
       call msg(msg_where,trim(ch),trim(db_ch(1)))
     else if (present(VAR)) then
       call msg(msg_where,VAR,trim(db_ch(1)))
     endif
   endif
   !
   ! DESCRIPTION
   !
   if (present(DESCRIPTION)) then
     DESCRIPTION=VAR//trim(db_ch(1))
   endif
   !
   io_rec_pos(ID)=io_rec_pos(ID)+len(db_ch(1))
 endif
 !
 contains
   !
   integer function i_verify(iv1,iv2)
     integer :: iv1(:),iv2(:),i1
     i_verify=0
     do i1=1,size(iv1)
       if (OP(i1)=="==") then 
         if (iv1(i1)/=iv2(i1)) i_verify=-1
       endif
       if (OP(i1)=="<") then 
         if (iv1(i1)>=iv2(i1)) i_verify=-1
       endif
       if (OP(i1)=="<=") then 
         if (iv1(i1)>iv2(i1)) i_verify=-1
       endif
       if (OP(i1)==">") then 
         if (iv1(i1)<=iv2(i1)) i_verify=-1
       endif
       if (OP(i1)==">=") then 
         if (iv1(i1)<iv2(i1)) i_verify=-1
       endif
       if (i_verify/=0) return
     enddo
   end function
   !
   integer function r_verify(rv1,rv2)
     integer  :: i1
     real(SP) :: rv1(:),rv2(:),eq_accuracy,acc_fraction
     r_verify   =0
     eq_accuracy=1.E-5_SP
     acc_fraction=0.001_SP
     !
     ! Define the accuracy of the comparison
     !
     do i1=1,size(rv1)
       if ( abs(rv1(i1))>0._SP) eq_accuracy=min( eq_accuracy , abs(rv1(i1))*acc_fraction )
       if ( abs(rv2(i1))>0._SP) eq_accuracy=min( eq_accuracy , abs(rv2(i1))*acc_fraction )
     enddo
     !
     do i1=1,size(rv1)
       if (OP(i1)=="==") then 
         if (abs(rv1(i1)-rv2(i1))>eq_accuracy) r_verify=-1
       endif
       if (OP(i1)=="<") then 
         if (rv1(i1)>=rv2(i1)) r_verify=-1
       endif
       if (OP(i1)=="<=") then 
         if (rv1(i1)>rv2(i1)) r_verify=-1
       endif
       if (OP(i1)==">") then 
         if (rv1(i1)<=rv2(i1)) r_verify=-1
       endif
       if (OP(i1)==">=") then 
         if (rv1(i1)<rv2(i1)) r_verify=-1
       endif
       if (r_verify/=0) return
     enddo
   end function
   !
end subroutine
